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1 Introduction 

A large number of effects related to the phenomenon of Bose-Einstein Condensa- 
tion (BEC) can be understood in terms of lowest order mean field theory, whereby 
the entire system is assumed to be condensed, with thermal and quantum fluc- 
tuations completely ignored [Ij. Such a treatment leads to the Gross-Pitaevskii 
Equation (GPE) [2] used extensively throughout this book. Although this the- 
ory works remarkably well for a broad range of experimental parameters, a more 
complete treatment is required for understanding various experiments, including 
experiments with solitons and vortices. Such treatments should include the dynam- 
ical coupling of the condensate to the thermal cloud, the effect of dimensionality, 
the role of quantum fluctuations, and should also describe the critical regime, 
including the process of condensate formation. 

The aim of this Chapter is to give a brief but insightful overview of various 
recent theories, which extend beyond the GPE (some of which may have already 
been mentioned in earlier chapters). As the author has been heavily involved in 
the development and implementation of such theories, the presentation will follow 
a somewhat personal perspective, following a logical progression from simple to 
more complex approaches, rather than a historical account. To keep the discussion 
brief, only the main notions and conclusions will be presented. While numerous 
appropriate key references are given, it is unfortunately not possible to mention 
all works that have contributed to the vast and growing field of finite temperature 
theories, and the author apologizes in advance for any such potential omissions. 

This Chapter generalizes the presentation of Chapter 1, by explicitly main- 
taining fluctuations around the condensate order parameter. While the theoretical 
arguments outlined here are generic, the emphasis is on approaches suitable for 
describing single weakly-interacting atomic Bose gases in harmonic traps. Inter- 
esting effects arising when condensates are trapped in double-well potentials and 
optical lattices, as well as the cases of spinor condensates, and atomic-molecular 
coupling, along with the modified or alternative theories needed to describe them, 
will not be covered here. 
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2 Microscopic Derivation of the Gross-Pitaevskii Equation 

For a sufficiently dilute ultracold atomic gas, where only binary interactions are 
relevant, the Hamiltonian of the system takes the form 

H = J dr#^(r,t)ffo<f'(r, ^) + \ j drdr'>P\r, t)4r'< {r' , t)V {r ~ r' )<F {r' , t)<F {r , t){l) 

where W^^'>{r,t) are the Bose field annihilation (creation) operators, ^(r — r') is 
the two-body interatomic potential, and Hq ~ —{h^/2m)V'^ + 14xt(r, i) is the 
'single-particle' operator in an external potential VJ;xt(r, t). In the usual discussion 
of BEC, the Bose field operator is typically split into two parts [3] 

!f'(t)(r,i) =!Z.(*)(r,i)cW(r,t)+V-'(t)(r,t), (2) 

respectively corresponding to the condensate, and non-condensate contributions. 
Here tf'*^*-' is the complex 'classical' amplitude of the condensate, while the operator 
(■•^t) is related to the annihilation (creation) of a particle from the condensate and 
accounts for fluctuations in the condensate phase. Initially, we restrict ourselves to 
relatively large three-dimensional (3D) systems, sufficiently far from the regime of 
critical fluctuations, where we can, to good approximation, ignore the effect of the 
condensate operator Q^^^ This approximation relies on the notion of spontaneous 
breaking of the U{1) symmetry, and leads to the appearance of a condensate mean 
field, !f'(r,t), often termed the 'condensate wavefunction', with all related fiuctu- 
ations contained in the 'non-condensate' operator ■0'(r, i). We will return to the 
role and physical interpretation of Q^\Y,t) at appropriate points in this Chapter. 

Symmetry- Breaking Picture: We expand the Bose field operator in terms of a com- 
plete basis set {0i(r)} via •f^^''') (r, t) — cl'P {t)4>i'^ (r), where brackets distinguish 
between creation and annihilation operators. We assume that each single-particle 
mode occupied by the condensate can be expressed in terms of a mean-field ampli- 
tude Zj-*^ = (a-^^), and a fluctuation operator c-^^ = a-^^ — (a-^^), obeying {c^p) = 0, 
where (• • •) denotes averaging over all states. The exact equation of motion for the 
condensate mean field amplitude of a general state n thus becomes [4j 



^2 ^ „ r 

ih—^ = hLJnZn + ^^{ni\V\jk) z*ZjZk + 2{c\cj)zk + z*{cjCk) + {c\cjCk) , (3) 

ijk 

where {ni\V\jk) = {l/2){{ni\V\jk) -\- {ni\V\kjy\ is the interaction matrix element 
with {ni\V\jk) = J J drdr' ip^{r)ilj*{r')V{r - r )ipj{r')ipk{r), and ?iw„ = {n\Ho\n). 
The first three terms in the sum of Eq. ([3]) depend on the condensate mean field, 
with the second term additionally depending on the non-condensate density, (cjcj ), 
and the third term on the average of two annihilation operators, (cjCk), termed the 
'anomalous average'. Importantly, however, Eq. (3) contains an additional term 
due to correlations of three fluctuation operators which does not depend directly on 
the condensate mean field. Although such 'triplet' terms (clcjCk) were ignored in 
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the majority of early treatments (see, e.g. [6l[71[8] and references therein), it is pre- 
cisely these contributions, which enable the transfer of atoms and energy between 
the condensate and the thermal cloud in kinetic theories [H [U [HI [IHl [HI HH US] ■ 

Introduction of an Effective Interaction: The usual 'naive' justification for the 
Gross-Pitaevskii Equation, is to ignore all non-condensate fluctuations, while si- 
multaneously approximating the interatomic potential by a delta function, via 
V{r — r') = gS{r — r'), where g a suitable effective interaction strength. While 
instructive, this is misleading, because V{r — r') corresponds to the actual inter- 
atomic potential, whereas g is the interaction strength for a complete effective 
binary s-wave scattering contact interaction in vacuum |14| , which generally pro- 
ceeds via virtual intermediate excited states, which may be occupied. This 'renor- 
malization' of the coupling constant was first discussed for the homogeneous Bose 
gas by Beliaev |15j and generalized to finite temperatures by Popov |16j; this 
problem was recently discussed by various approaches by Stoof et al. [TTl [18], 
Proukakis et al. [HI [201 [21] , Shi and Griffin and Morgan [23]. We follow the 
notation of Proukakis et al. [5] and investigate the condensate evolution at zero 
temperature, where there is a negligible non-condensate population, such that 
(cjcj) « {c\cjCk) ~ 0. Apart from pure mean field effects, z*ZjZk, the equation for 
the mean amplitude z„ still contains a contribution from the anomalous average. 
The latter can be separated via (cjCk) — {cjCk)Q+6{cjCk) into a slowly-varying con- 
tribution over low-lying states, {cjCk)o, and a contribution over high-lying states, 
S{cjCk), which varies within a coUisional timescale. By using the corresponding 
equation of motion for S{cjCk) to adiabatically eliminate this quantity from the 
zero-temperature limit of Eq. ([3]), the latter takes the form [4l [T9] 



where T represents an effective two-body interaction term obeying a Lippmann- 
Schwinger relation [|4; expressing the complete repeated scattering amplitude for 
two particles in an empty trap. This procedure corresponds to a diagrammatic 
resummation of a selected class of diagrams |3J [T71 [S^ . This effective interaction 
is strictly only valid when used to describe the coupling between low momen- 
tum states [H [131 [531 121] j ^ind cannot be used as a real potential acting between 
the atoms. Ignoring all energy and momentum dependence, this 'renormalized' 
interaction can now be approximated by an effective contact (pseudo)potential 
T^^(r, r') = gS{r — r'), where g — Anh^a/m, and a is the s-wave scattering length 
[14[ 125 ]. The implicit assumption of a constant, momentum- independent, scatter- 
ing amplitude must be used cautiously, and is valid only up to some maximum 
momentum /cmax oc 1/a |16j . The anomalous contribution {cjCk)o accounts for 
many-body effects which are typically neglected in dilute weakly-interacting gases 
(see Sec. [3]). Re-expressing the mean- value equation in terms of the condensate 
mean field !f'(r,t), we obtain the Gross-Pitaevskii Equation at T = 0, 




huJnZn + ^{ni\f'^^\jk)z* [ZjZk + (CjCfe)o] , 



(4) 



ijk 



ih-^{v,t) = \H^+g\^{v,t)\'] tf'(r,i) =i/GP-f(r,t) 



(5) 
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Number-Conserving Treatment: The above treatment assumes the existence of a 
mean field, and thus rehes on spontaneous symmetry breaking. While convenient, 
this notion is somewhat ill-defined for a finite trapped gas with a fixed atom 
number iV, since the condensate mean field corresponds to the ensemble average 
over number states |iV), with {N\^{r,t)\N) = <I'{r,t){N\N - 1) = 0. The usual 
justification is that for a system with a large number of particles, N ±1 N with 
accuracy (1/A^), essentially assuming that the condensate operators commute [26]. 
Such arguments are valid in the grand canonical ensemble, where the atom number 
is not fixed, and one minimizes the shifted hamiltonian {H — fJ,N), where /i the 
chemical potential and N the number operator. Also, no difficulties arise if, by 
analogy to laser theory, we treat the condensate as a coherent state |a), i.e. a 
superposition of number states for which {a\a) 7^ by construction. 

Gardiner argued that the apparent non-conservation of particle numbers arises 
from an inappropriate definition of the quasi-particle operators. This led to the 
development of a C/(l)-symmetry-preserving approach, relying on the separation 
of the 'particle' and 'phonon' concepts [2I]. Following Gardiner's notation, 'f'(r) = 
A[^(r)-|-(l/-\/]V)x(r)] where ^(r) is the condensate wavefunction obeying the usual 
GPE, x(r) is a phonon field operator, and A is a generalized annihilation operator 
ensuring total particle number conservation. This approach sets the correct basis 
for further treatment [28], with an equivalent approach formulated by Castin and 
Dum [5S] and extended to finite temperatures by Morgan [53]. 

The remaining sections generalize the above treatment to finite temperatures. 
Crudely speaking, one can divide such generalized treatments into those that rely 
on a suitably truncated hierarchy of coupled equations of motion for generalized 
mean fields (Secs.[3]&|4|), and those in which the lowest modes of the system, which 
include both condensate and quasi-particlcs, are described in a unified manner 
by appropriate probability distribution functions (Sec. [5]), giving rise to various 
stochastic approaches (Sec. |6|). Sec. [7| briefly discusses the importance of fluctua- 
tions in the phase of the condensate, relevant in low dimensional systems. 

3 Generalized Mean Fields: Static Thermal Cloud 

From Eq. (3), it is clear that the condensate mean field !f'(r, t) is coupled to 'higher 
order' mean fields (normal and anomalous averages). Their corresponding equa- 
tions of motion are, in turn, coupled to yet higher order correlations of the non- 
condensate operators, and so forth. A closed set of equations requires a consistent 
set of decoupling approximations, thereby uniquely defining a set of generalized 
mean fields which should accurately determine the static and dynamic properties 
of the system }30j . Equations of this form are well-known in the literature and can 
be derived by variational [17[I31| . or diagrammatic ^22j methods. A description in 
terms of the three lowest order averages {'P), {ip'^'ip') and {ip'ip') is known as the 
Hartree-Fock-Bogoliubov (HFB) theory (see below) [32] , 
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Mean Field Approximation: To first approximation, one can ignore particle ex- 
change between condensate and thermal cloud (i.e. the triplets), and allow only 
mean-field coupling between them. Such a theory arises upon reducing the hamil- 
tonian of Eq. ([T]) to quadratic form via the 'mean field' approximations [H [S] |23] 
i^'^'-ip' « 2{'ip'^')^' + and « 4(V;'tV'')^'^-i/'' + (?/''1'?A'1')V''V'' + 

(V'V')V''^^'^ - 2(V''^')(^'^^'> - the latter justified from Wick's 

theorem |14) . This leads to an effective quadratic hamiltonian in which single- 
particle energies are shifted by the added mean field potentials. Such a hamilto- 
nian can be diagonalized by the Bogoliubov transformation ?A'(r) = 0i(r)ci = 
^j[uj{v)j3j — i;*(r)/3j], mapping single-particle non-condensate operators Ci into 
mixtures of 'quasi-particle' operators j3j and /3j, This leads to a generalized GPE 

ih^^^{r,t)^ [i/o + ff [|'^'(r,t)P + 2n'(r,i)]] if(r,t)-f mo(r,t)f*(r,i), (6) 

which, compared to the GPE of Eq. ([5]) contains additional contributions from the 
mean field coupling of the condensate ^ to the non-condensate n' (with the factor 
of two arising due to direct and exchange contributions [14j ) and the anomalous 
average mp. The above equation is coupled to finite temperature 'Bogoliubov-de 
Gennes' equations for the excitation amplitudes Uj{r), Vj{r). Such a closed set of 
equations (see Eqs. Q-®) is known as the HFB theory [HI [HISIl 132 ■ 

Many-Body Effects: We argued that the two-body effective interaction, T^^, arises 
via adiabatic elimination of rapidly-varying anomalous averages. Such an effec- 
tive interaction describes the coUisional amplitude in vacuum, whereas in trapped 
gases collisions take place within a medium of other condensed and excited par- 
ticles. This leads to important modifications in the description of atom-atom in- 
teractions [171 im EQI [22] which can again be incorporated into a more general 
effective interaction, the many-body T-matrix, T^^. This can be approximated 
(in the zero-energy, zero-momentum limit) by a contact potential with a position- 
dependent amplitude .gcff(i") = 5[1 + "^o(r)/'0^(r)] [HI [211 [M] ; here mo(r) denotes 
the 'regularized' position representation of the anomalous average over low-lying 
states (see below). Since mo(r) is of the same order of magnitude as the non- 
condensate density n'(r), its effect should, in general, not be ignored [33] . 

HFB and Related Approaches: To first approximation, the non-condensate com- 
ponent is treated as static, yielding the static HFB equations, which, written in a 
slightly more general notation to facilitate subsequent discussion, take the form 

^o-/i + 5co„(r)|V^(r)|2 -f 2gexc(r)n'(r)l ^(r) = (7) 



-M* -L* J \ Vj J ^ \ Vj ' ^ ' 



where L{r) = Hq ~ ^l -\- 25con(r)|'(A(r)h + 25cxc(r)ri'(r) and M = 5con(r)V' (r 
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This general form distinguishes between the effective interaction for a collision 
between two condensate atoms, (7con(r), from the corresponding one for collisions 
involving a condensate and a non-condensate atom, goxc(r). Solving Eqs. Q-® 
yields the 'Bogoliubov functions' Uj{r) and Vj(r), in terms of which we can directly 
obtain the time-independent normal n'{r) = + bj(i")P)^BE(ej) + 

|fj(r)p and anomalous mo(r) — J2j ^j(i")''^j'(r)[2A^BE(£j) + 1] — TOuv(r) averages. 
Here A^BE(£j) = [exp {(3hujj) — 1] ^ is the Bosc-Einstcin distribution with /3 = 
l/ksT. The term muv(r) is required to ensure the regularization of the anomalous 
average by subtracting off the ultraviolet divergences arising in the homogeneous 
limit due to the use of a momentum-independent contact interaction. 

The standard HFB limit is based on (/con(r) = goxc(r) = g everywhere except in 
Eq. ([7]) where (7exc(r) — g still, but 5con(r) — 3cff(r). This inconsistent treatment of 
effective interactions leads to the appearance of a gap in the homogeneous single- 
particle excitation spectrum at zero momentum, in direct contrast to the Goldstone 
theorem which guarantees a gapless spectrum in the presence of U{1) symmetry- 
breaking [35], and the Hugcnholtz-Pines theorem [36]. As a result, the full HFB 
theory should not be used to calculate excitation frequencies. One common way 
to overcome this problem, is to completely ignore the additional influence of the 
medium on collisions of trapped atoms, i.e. to set mo(r) = [H [22]. This led to 
the first quantitative predictions of excitation frequencies of trapped Bose gases at 
low temperatures 371 , in agreement with experiments. The reduced gapless theory 
with mo(i") = is often termed 'HFB-Popov' approximation (but see [33].) 

In order to obtain a consistent theory for excitation frequencies, but addition- 
ally include the modification imposed by the medium on the collisional properties 
(i.e. mo(r 7^ 0)), one may use the above generalized form of the theory, with the 
effective interaction gcon(r) = ffcff(r) for the scattering of two atoms in the con- 
densate, and the mean-field coupling between the condensate and non-condensate, 
5cxc(r), approximated either by or by goff(r) consistently throughout Eqs. ([7])- 
([8]) . These correspond to different approximate treatments of the relative momenta 
of colliding atoms, and give rise to two different generalized 'gapless' theories, 
GHFBl and GHFB2 [H [Mj- (see also [3H]). These theories have been applied 
to excitation frequencies, vortices and the study of coherence of 2D condensates 
[34] [39] . Alternative approaches to the HFB gap problem appear in [33] i40] ■ 

4 Generalized Mean Fields: Dynamic Thermal Cloud 

So far, the non-condensate has been treated as static. Such approaches cannot 
account for common dissipative process in the condensate, e.g. Landau damping 
(coalescence of two excitations into a single one of higher energy) and Beliaev 
damping (breaking up of an excitation into two excitations of lower energies) , which 
require a coupled dynamical description of both condensate and non-condensate 
[T] . A complete theoretical picture must also describe collisions between condensate 
and non-condensate atoms, as well as between two non-condensate atoms. 



Beyond Gross-Pitaevskii Mean Field Theory 



7 



4.1 Time-Dependent Hartree-Fock-Bogoliubov 

This section outlines the work of the author [3] as finalized in [5^. The exact 
equation of motion for the non-condensate {c\cj) depends on averages of up to four 
fluctuation operators Ci — Ui — Zi. Using second order perturbation theory, and 
considering for simplicity a single condensate mean field amplitude, zq, coupled to 
numerous excited levels of populations rii , one obtains coupled equations governing 
their evolution. The equation for dui/dt includes, among other terms, the following 
contribution accounting for collisions between two non-condensate atoms 

[j2]'^\^rfrni\ [(j^i + 1) (?^m + - "iMm + 1) (»^s + 1 )] ^(t^rsmi ) (9) 

^ ^ rms 



where ujrsmi ~ i'-^r + t^s — (^rn ~ ^i) ■ We also obtain the following additional con- 
tribution to the thermal cloud dynamics drii / dt due to the presence of a condensate 



^) \zo\'Y.\T??o.\\in^ + l)r 



rUs - ni{nr + l){ns + 1)] 5{uJrsOi)- (10) 



This term is proportional to |zoP and describes the change in the non-condensate 
population of level i induced by the transfer of a single atom to, or from the 
condensate (level 0). One must also consider the corresponding change in the 
condensate population associated with such processes. Adiabatic elimination of 
the triplet correlation {c\cjCk) in Eq. ([3]) leads to [4] 



zo (11) 



+ i (-^j ^0 ^ l^i^s^if li^i + l)«r?^s - f^i(?^r + l)("-s + 1)] 5 (uJrsO'i.) 
irs 

which clearly contains both energy shifts due to mean fields and damping. The 
equations given above contain the essential key elements for a consistent kinetic 
theory for a partially-Bose-condensed gas. These reveal that while the scattering of 
a particle into an excited state i is 'bosonically enhanced' by the factor (rii + l), the 
corresponding scattering into the condensate does not feature spontaneous growth 
(i.e. one simply obtains |zoP and zq in Eqs. pil)) and (fTTj) respectively). 

The discussion thus far relied on a perturbative expansion beyond a single- 
particle basis of bare (undressed) particles, and has led to dressing of the energies 
from the mean fields of the condensate, |^oPj and the non-condensate, n^. Such ba- 
sis choice is, however, not very suitable for typical systems containing a condensate 
spanning numerous trap eigenstates. While one could transform the above equa- 
tions to a quasi-particle basis, this leads to rather involved expressions [41] |42] . 
An alternative, but equivalent, approach is to first incorporate all relevant HFB 
mean field effects into the basis, i.e. shift from bare single-particle energies TiLUn = 
{n\Ho\n) to dressed quasi-particle energies hui^^^ = {n\HiiPB\n) where Hhfb is 
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the reduced quadratic hamiltonian with mean field effects included [14j, and sub- 
sequently perform second order perturbation theory beyond this dressed basis. In 
this case, the perturbative terms arise from the difference between the exact ex- 
pressions of products of fluctuation operators, and their corresponding HFB mean 
field decompositions. For example, one part of the perturbative hamiltonian arises 
from the 'beyond-mean-field' quantity {iP''<iP'^'i^')-(2{'iP''fi}')^ + (?/>'tV^'t) (?/,'^')^ 

as suggested by Morgan [23] and Zaremba, Nikuni and Grifhn [43 . Explicitly 
performing this (see for details) yields a set of coupled time-dependent non- 
markovian equations for the mean field amplitudes z, (c^c) and (cc), corresponding 
to n' and rfiQ. These equations describe all atomic collisions, including not only 
collisions within the condensate, but also particle exchange between condensate 
and non-condensate and collisions of two thermal atoms. This theory was first 
derived by Walser et al. [12l|44], using slightly different, but entirely equivalent, 
arguments (based on separation of timescales for the evolution of the relevant 
physical parameters), and subsequently obtained by the author 9]. The reader is 
referred to [Hllllllll] for details, implementation, and relation to other treatments. 

The perturbation theory presented here can also be derived, perhaps somewhat 
more systematically, from the method of non-commutative cumulants, in which the 
cumulants capture the essential correlations in the system by subtracting factoriz- 
able contributions from the correlation functions, as mentioned above, with such 
an approach already successfully implcmcted in diverse contexts [45j . 

Gapless Second- Order Perturbation Theory: To monitor situations close to equi- 
librium, Giorgini and Rusch and Burnett [46| performed a linear response 
treatment around static values of the HFB parameters z, (c'c), and (cc). This was 
discussed in detail in a series of papers by Morgan ^23 j , who additionally reconciled 
the rigorous introduction of an effective interaction with the number-conserving 
formalism. Morgan's approach, based on a perturbative treatment beyond the 
HFB basis, accounts for the indirect excitation of the condensate via excitation of 
the thermal cloud, as needed to interpret the excitation frequencies and damping 
rates at JILA [47] (see also [48| [49]). The resulting form of the Bogoliubov-de 
Gennes equations for the excitations are similar to those of Eqs. (O-®, but with 
additional, or modified terms to ensure orthogonality between condensate and non- 
condensate and number conservation. Such a theory consistently includes damping 
effects, population changes due to quasi-particle collisions and finite size effects. 

4.2 Theory of Zaremba-Nikuni-GrifRn 

Parallel to the above developments, Zaremba, Nikuni and Griffin (ZNG) presented 
the first computationally implementable generalized mean field theory describing 
the coupled evolution of condensate and thermal cloud [43] . The starting point is 
the same as the time-dependent HFB approach. While the subsequent treatment 
ignores the anomalous average (i.e. many-body effects included in [HI [H]), this 
theory made a significant contribution by implementing local energy and momen- 
tum conservation [30j . The resulting kinetic equations are valid in both mean-field 
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and hydrodynamic regimes. In this theory, the condensate evolves according to 
d 



Ho+g[\^{v,t)\^ + 2n\v,t)\-iR{v,t) ^{v,t) (12) 



with the coupled non-condensate density n'{r,t) described in terms of a Wigner 
phase-space representation /(r, p, i) via n'{r,t) = J dp/(27r/?i)'^/(r, p, i), with / 
obeying a generalized 'Quantum Boltzmann' kinetic equation [501 

-/(r, p, t) + ^ • V/(r, p, t) - VX^eff (r, t) ■ Vp/(r, p, t) = Ci2 [/] + C22 [/] (13) 
ot m 

with an effective mean-field potential Vefi(r, t) — Vcxt(r, t)+2g [\^{r, t)p + n'(r, t)] . 
This equation contains two coUisional terms: C22 describing collisions between 
non-condensate atoms, and C12 describing the transfer of an atom from the non- 
condensate into the condensate, and its inverse process. These are given by 

y (2^ y (2^ y 12^ ^'^'^ (p + p. - p., - P4) 

5{e + e2~e3~ £4) [(1 + /)(! + /2)/3/4 - 7/2(1 + /3)(1 + A)] (14) 



S {ec + S2-e3- £4) (27r?i)^ [S{p - P2) - S{p - ps) - 6{p ~ P4)] 

X [(1 + /2)/3/4 - /2(1 + /3)(1 + /4)] , (15) 

where fi = f{ri,pi,t), Si = pf/2m + T4fr(r), and Ec = + mv1/2 is the local 
condensate energy. The condensate growth term, —iR{r, t) in Eq. p2|l arises from 
'triplet' correlations and is defined by i?(r, i) = h/2\<I'\'^ J dp/{2Trn)^Ci2[f]- 

The original numerical implementation of this theory jllj successfully ad- 
dressed the issue of condensate formation by relying on the ergodic approximation 
[50] , which assumes that equilibration is rapid for atoms with similar energies, thus 
enabling the phase-space variable fi = f(ri,pi,t) to be expressed only in terms 
of an energy variable E{ri,pi,t). These calculations were subsequently general- 
ized [m El] to strongly non-equilibrium regimes by representing the phase-space 
density by a collection of N discrete 'test particles', with collisions between them 
handled via an appropriate Monte Carlo sampling scheme. Despite ignoring many- 
body effects, this theory yields remarkable agreement with various experiments, 
e.g. scissor's mode, quadrupole excitations and transverse breathing modes of elon- 
gated condensates jSlj. Hydrodynamic collective modes and vortex nucleation at 
finite temperatures were also discussed [SH [S3] . The ZNG theory also explained 
the observed dark soliton decay in the Hannover experiment JS4J, and ascertained 
the required conditions for the observation of dark soliton oscillations in elongated 
harmonically-confined BECs [55] . The ZNG theory arises as a limiting case of the 
approaches of Stoof and Gardiner, ZoUer and co-workers discussed next. 

Damping of excitations was also addressed by a different approach by Fedichev, 
Shlyapnikov and collaborators who studied stability of vortices and dark solitons 
by considering the scattering of excitations from such macroscopically excited 
states [57j , with this formalism applied to bright solitons by Sinha et al. [58] . 
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5 Kinetic Theories based on Probability Distribution 
Functions 

All above theories have the common feature of treating the condensate as a sep- 
arate entity from the thermal cloud, with the former typically identified as the 
mean value of the Bose field operator {^{r, t)). Such perturbative theories cannot 
predict phenomena in the critical regime, e.g. onset of condensation, or the shift in 
Tc- An alternative 'class' of approaches was thus developed, in which one treats the 
low-lying part of the spectrum, including both the condensate and the low-lying 
excitations as one entity, with the remaining high-lying modes treated separately, 
by means of a Quantum Boltzmann Equation. 



5.1 Stoof s Non-Equilibrium Theory 

Stoof developed a unified finite temperature quantum kinetic theory by investigat- 
ing the evolution of the entire system by means of a Wigner probability distribution 
P[4>* t] [T8l[56]. At time t, this distribution gives the probability of the system to 
be in a coherent state |(/'(r); t) = exp{/ (ir0(r)tf'T(r, i)}|0), where |0) is the vacuum 
state. Using suitably normalized coherent states, P[(f>*,(t); t] = Tr [p(to)|0; t) (0; t\] , 
where /3(to) the initial density matrix of the system prior to the onset of con- 
densation, the determination of the probability function reduces to the calcula- 
tion of a functional integral containing the probabilities | (0; t|(/)o; to) This can 
be evaluated as a 'path integral' over all field evolutions from t to to and back 
to t, leading to the Keldysh non-equilibrium formalism. Such an approach, en- 
ables a systematic consideration of the fiuctuations around the mean field, with 
all correlation functions obtained from moments of the distribution P[(j)* ^(p^t]. 
The derivation of this theory is given in 56J , with a brief insightful account pre- 
sented in [59l [60l [61] . By working with a quantized quadratic effective action, in 
which interactions are 'renormalized' to the many-body T-matrix, one obtains a 
Fokker-Planck equation for the full probability distribution function P[0*,(/);t]. 
One then 'decouples' condensate and non-condensate modes, via the substitution 
P [(j)* ,cj)] t] = Fo [^*, t] Pi [(j)'* , cj)' ; t\. This gives rise to a coupled set of equations: 
firstly, a Fokker-Planck equation for the temporal evolution of the probability dis- 
tribution of the condensate wavefunction in the presence of a thermal cloud [56] 

z^^Po [1>\1>-A = - J drj^^ [Ho - m - iR[v,t)+g\<l>[v)\^'\ <P{r)Pa [<P*,<P;t] 

+ J drj^ [Ho~Kt)+iRir,t)+g\<Pirt] <P* {r)Po ,<P;t] 

- \ J -^- ^.^(r)l*(r) '^*(-)"^"(-'^)^°['^*''^^^]- 

Here zi?(r, t) corresponds to particle exchange between the condensate and non- 
condensate, while the Keldysh self-energy {r, t) (see Eq. ((T9)) ) expresses thermal 
fluctuations due to incoherent collisions between condensate and non-condensate 
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atoms. The above equation contains both 'streaming' and 'diffusion' terms; the 
absence of third-order derivatives (present in [5H arises from the con- 

sistent ehmination of selected terms in the effective action, required to avoid 
double-counting in the many-body approximation. The non-condensate dynam- 
ics is parametrized by introducing a Wigner distribution via TV (r, p, i) -I- 1/2 = 
/dr'e"'P '"'((/)'(r + r72)0'*(r-r72))(i) with iV (r, p, i) obeying a Quantum Boltz- 
mann Equation. The formulation of this theory ensures that the fluctuation- 
dissipation theorem relating R{r,t) and S^(r,t) is satisfied at equilibrium, thus 
guaranteeing that the condensate relaxes to the correct equilibrium distribution. 
Although qualitative predictions regarding condensate formation can be made, the 
solution of this reduced Fokker-Planck equation under general conditions is compli- 
cated, because terms appearing within this equation depend implicitly on <P'^*^ via 
their dependence on energy. One therefore typically resorts to approximate solu- 
tions, by means of the stochastic GPE discussed in Sec.[6l Related non-equilibrium 
quantum field theoretic approaches have been recently discussed in |64j . 

5.2 The Gardiner-Zoller Quantum Kinetic Master Equation 

In a series of papers [65], Gardiner, ZoUer and co-workers used techniques from 
quantum optics to derive a quantum kinetic master equation for the condensate, 
coupled to a time-dependent thermal cloud. This treatment is based on splitting 
the physical system into two 'energy bands', each of which has its own evolution, 
but with the two bands additionally exchanging both particles and energy: (i) 
The 'Condensate Band', Rc, is a band of low-lying states which typically include 
the condensate, and those low-lying non-condensate modes that are heavily af- 
fected by the presence of the condensate. It is formulated within the Bogoliubov 
number-conserving approximation, with the basis expansion performed in terms 
of many-body eigenfunctions. This band is fully described by the total number of 
atoms and the quantum state of the quasi-particles within Rc ■ While the atoms are 
conserved, the quasi-particles are mixtures of phonon states, and their number may 
change by the absorption, or creation, of a quasi-particle in a single coUisional pro- 
cess, (ii) The 'Non-Condensate' Band, Rnc- Atoms with energies higher than some 
appropriate cut-off energy, En, are described in terms of a wavelet expansion, with 
particles binned into discrete shells in phase space, on the basis of their position and 
momentum. This is often treated as being in thermal equilibrium, with a local par- 
ticle phase-space density /(r^, p^, t), temperature T, and chemical potential ii. The 
master equation for the condensate band, including both condensate and quasi- 
particles, is thus obtained [65] . following the standard methodology of quantum 
optics _66J. A simplified version gives rise to a rate equation for the mean number, 
Nc, of atoms in the condensate [28| dNc/dt = 2W+{Na)[{l - exp(-/3Z\/i)7Vc + 1], 
where W'^{Nc) ~ ['im{akT)'^ /nh^] e^^^ denotes a growth rate depending on the 
system scattering length a, temperature T, and chemical potentials of the conden- 
sate and non-condensate, with Z\/x the difference in chemical potential between 
the two subsystems. This approach led to the first quantitative predictions of con- 
densate growth [28] , and to good agreement with experiments [6ll [68] ■ Subsequent 



12 N.P. Proukakis 



improved work also took into account the evolution of occupations of lower trap 
levels via the ergodic quantum Boltzmann equation. 

6 Stochastic Approaches to Condensate Dynamics 
6.1 Classical Field Methods 

Although the preceeding discussion indicates severe limitations to the validity of 
the GPE at finite temperatures, remarkably the GPE can actually be used to model 
the evolution of any system described by the Hamiltonian of Eq. (1) provided that 
such a system primarily behaves in a classical manner. Such classical evolution was 
first studied in the context of condensation by Kagan, Svistunov and Shlyapnikov 
[69j . This approximation was also used to study an ideal superfluid approaching 
equilibrium in [70], with a qualitative 2D study of evaporative cooling performed 
in [71] , and large scale homogeneous simulations performed in [72] . The use of the 
GPE to study a weakly-interacting trapped gas stems from the realization that 
when all relevant modes of the system are highly occupied {uk S> 1), the classical 
fluctuations of this field are much larger than quantum fluctuations. Such modes 
can then be fairly accurately represented by a coherent wavefunction, in analogy 
to the classical description of highly-occupied laser modes. 

The Projected Gross- Pitaevskii Equation: A classical description can only work for 
modes up to a certain energy cut-off, since the condition of high occupation will 
break down for high-lying modes. Such a cut-off can be implemented through the 
use of a projector, P, onto highly-populated modes, giving rise to 

i?i^tf'(r, t) = Ho^iv, t)+gP{ |tf'(r, t)\^ 9{r, <)} + /(t}). (17) 

The term /(?)) appearing above is a rather complicated function of the 'fluctuation 
operator' = (1 — P)!f'(r) [TSIIT^, describing the coupHng of the classical region to 
a heat bath, whose modes are not highly populated and are orthogonal to the con- 
densate. Ignoring /(?)) leads to the 'Projected Gross-Pitaevskii Equation' (PGPE) 
introduced by Davis, Morgan and Burnett [73]. This is similar to the usual GPE, 
Eq. ([S]), apart from the presence of the projection operator, P, which ensures that 
all computed quantities remain within the classical region. The projector used in 
this theory is diagonal in the single-particle basis of the hamiltonian and is de- 
fined by P{P(r, t)} = J2nec <^"(^) / '^^rc/)* (r')P(r', t) where (^„(r) denotes the 
eigenf unction, with n restricted within the coherent (classical) region C. Although 
the PGPE contains no damping terms, since the quantity /(f)) has been dropped, 
the presence of elastic collisions ensures that a highly non-equilibrium initial con- 
figuration relaxes to an equilibrium distribution of given energy, without being 
particularly sensitive to the initial conditions. Such a distribution can be assigned 
a temperature, upon noting that the mean occupation will actually be given by the 
classical limit n{ei) = [P (e^ — /z)] ^ of the Bose- Einstein distribution. Importantly, 
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the field 'P{r, t) obtained from the PGPE represents the quantum field of many 
low- lying modes, rather than a single condensate mode. This theory has been used 
to investigate the evolution of vorticity in homogeneous systems and generalized 
to trapped systems to discuss the process of evaporative cooling, the shift in Tc 
and the study of spontaneous vortex-antivortex pair production in quasi-2D gases 
[74j . Since the PGPE is non-perturbative, it can be applied even at the critical 
region. Related classical field approaches studying the equilibration of a suitably 
randomized initial condition of specified energy by means of the GPE appear in 



The Truncated Wigner Approximation: The idea of evolving randomized initial 
conditions via a classical field equation, such as the GPE, can be rigorously jus- 
tified by consideration of the Wigner quasi-distribution function. In general, the 
equation of motion for the Wigner distribution of a gas of ultracold atoms con- 
tains both first and third order derivatives with respect to the atomic field, and the 
Truncated Wigner Approximation (TWA) amounts to ignoring the latter higher- 
order derivatives. This is justified provided that the number of particles in the 
system is much larger than the number of accessible modes, such that all modes 
have an occupation larger than unity. Thus, the main idea behind this approach 
is to create a set of random classical fields which accurately sample the Wigner 
distribution function of the density operator of the system, and evolve these via 
the Gross-Pitaevskii equation. Despite the deterministic nature of this equation, 
effects of quantum noise are actually maintained in an approximate manner in 
the initial state, in the form of mode amplitude fluctuations |62j . which should be 
appropriately sampled. The first numerical implementation of the TWA was per- 
formed in ,63J, with appropriate sampling techniques analysed in |76j . The TWA 
was recently used to study colliding condensates condensate reflection from 
a steep barrier, three-body recombination processes, and collapsing condensates 
|77| . Polkovnikov performed a systematic perturbation theory in quantum fluctu- 
ations around the classical system evolution, obtaining the GPE to lowest order 
and the TWA to next order W8\. 




In addition to the (truncated) Wigner representation, a quantum system de- 
scribed by a master, or Fokker-Planck, equation, can also be studied within the 
context of the positive-P representation [M , which gives exact results provided 
that the ensemble averages converge, thus restricting both system size and evo- 
lution timescales. The evaporative cooling dynamics were studied by Drummond 
and Corney J79] , with significant progress made recently f8(T . This discussion omits 
various related techniques, e.g. Monte-Carlo methods. 

In our treatment so far, the system relaxes to some equilibrium state, due 
to noise implemented in the initial conditions only. However, the only way to 
guarantee that the system will relax to the correct equilibrium is to allow the 
condensate to be coupled to a thermal cloud, which provides in an intuitive manner 
the necessary irreversibility. One possibility to extend beyond the PGPE or the 
TWA is to treat the high energy modes as a thermal reservoir, as discussed below. 
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6.2 The Stochastic Gross-Pitaevskii Equation 

The Fokker-Planck equation derived by Stoof [56 , Eq. (|16|) . can be mapped directly 
onto a Langevin equation [66j for the order parameter of the system ^(r, t), with 
muhiphcative noise |61| . This problem can be simplified considerably by assuming 
that the thermal cloud is sufficiently close to equilibrium, that it can be described 
by a Bose distribution with a chemical potential fi and temperature T. In this 
regime, the thermal cloud acts as a heat bath to the condensate, and the Langevin 
equation for the condensate takes the form [5T] 

i?i^^>(r, t) = [(Fgp -f^) - ii?(r, t)] <^(r, t) + ry(r, t) (18) 

where Hgp = Hq + g\<P{r,t)\'^ is the usual 'Gross-Pitaevskii' hamiltonian. This 
Stochastic GPE (SGPE) describes both mean field effects of the condensate and 
low-lying excited states, and fluctuations about the mean field, thus providing valu- 
able information on both diagonal and ofF-diagonal elements of the one-particle 
density matrix. In the 'classical' approximation, where N{ei) = [P{£i — fi)] ^, 
the dissipative term obeys iR{r,t) = —{f]/A)'hX!^{r,t)HGP- The contribution 
ri{r,t) denotes a 'noise term', with Gaussian correlations {rj* {r,t)ri{r' ,t')) — 
i{n^ /2)E^{r, t)S(t - t')S{r - r'). The equivalence between Eq. and Eq. p6|) 
comes about because averaging the product <?*(r, t)<?(r, i) over the different re- 
alizations of the noise ry(r, t) is by construction equivalent to averaging over the 
Wigner distribution P[<?*,<?;t]. The Keldysh self-energy E^{r,t) describing the 
strength of the thermal fluctuations due to incoherent collisions between conden- 
sate and non-condensate atoms is given by 

4:n\ 2 f dp2 f dp3 f dp4 3 



6 {Ec + 62-63- £4) X [(1 + N2)N3N4 + N2{1 + N3){1 + N4)] (19) 

where Ni = N{ei) the Bose distribution for the eliminated part of the gas, with en- 
ergy Ei = (pf -I- T4xt(r) + 25|^P), and Sc corresponds to the average local chemical 
potential. A simplified form of this SGPE, assuming time-independent occupation 
numbers in the non-condensate, and thus a time- independent self-energy (r) 
has been used to study [51] the reversible condensate formation observed experi- 
mentally when cycling through the phase transition |82| , quasi-condensate growth 
on an atom chip, fluctuations of one-dimensional Bose gases, and the growth of 
coherence of an atom laser [53]. This theory is also amenable to variational cal- 
culations, and such a technique has been used to discuss coUisional frequencies 
and damping rates of collective excitations, growth-collapse cycles in attractive 
condensates [61) , and finite temperature dynamics of a single vortex [84] . 

Additionally ignoring the 'noise term' 77 in Eq. (|18p . in the 'classical' approxi- 
mation discussed here, the condensate evolution can be cast in the form 

z?i^<?'(r, t) = (1 - Z7) (Hgp - <?(r, t), (20) 
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with a temperature- and position-dependent damping rate 7 = z(/?/4)?iZ'^(r). 
Phenomenological damping of this form was originaUy proposed by Pitaevskii |85| , 
and first implemented to trapped Bose gases with a constant, position- independent 
rate 7 by Choi et al. [86j . A related phcnomcnologically-damped equation with the 
factor (1 — 17) appearing on the left hand side of the equation has been used in 
diverse studies, including vortex lattice growth and dark soliton decay [87] . 

A similar SGPE has been recently derived by Gardiner et al. [88] in terms of the 
Wigner function representation. This combines the ideas of the quantum kinetic 
theory of Gardiner, ZoUer and co-workers [65] with those of the finite temperature 
GPE [13]. The system is formally split into low- and high- lying modes by means 
of the projectors P and (1 — P), with the high-lying modes treated as a heat 
bath. This leads to a master equation for the interaction of a condensate with a 
fixed bath of non-condensed atoms, which can be mapped onto a Fokker-Planck 
equation. Subtle differences arise between the approaches of Stoof and Gardiner, in 
the effective interactions (many-body vs. two-body), the elimination of high energy 
modes, and the introduction and type of noise terms. When projection operators 
and noise contributions are ignored, this equation leads precisely to Eq. (|20p . which 
has been used to study rotating condensates and vortex lattice formation [89] . 

7 The Role of System Dimensionality 

In a homogeneous system, BEG arises only in 3D (or 2D at T = 0). Recent ex- 
perimental progress has enabled the achievement of weakly-interacting quasi-2D 
(quasi-lD) geometries, in which motion in the transverse direction(s) is 'frozen out' 
by tight harmonic confinement, whose typical level spacing ?iuj± exceeds both ther- 
mal and interaction energies. Although condensation is approximately regained in 
such reduced dimensionality systems at sufficiently low temperatures, a partially- 
condensed system, in general, suffers from phase fluctuations which limit the co- 
herence to a region smaller than the 'condensate' spatial extent, as experimentally 
observed ^U\; the system is then said to contain a 'quasi-condensate' [T^. While 
the stochastic theories presented earlier can be used to describe such systems [S^ , 
any attempt to perform a mean field analysis of this regime requires consideration 
of the 'condensate operator' ( of Eq. ([2]) [91]. In general, the Bose field operator 
should be expressed as 'P{r) — ^no(r)e**('") -I- '0'(r) [92^, where, uq is the quasi- 
condensate density, and x(r) an operator accounting for phase fluctuations. Clearly 
no must vanish in the thermodynamic limit in ID and 2D (T 7^ 0) and yield the 
condensate density in 3D, which occurs when the condensate density is identi- 
fied as Lim {noexp{— (l/2)[x(r) — x(r')]^}} when (r — r') — + 00. This approach 
provides ab initio predictions for the equilibrium properties of weakly-interacting 
trapped gases in all dimensions, in agreement with the stochastic GPE [83] . 

As soon as phase fluctuations come into play, the GPE fails to predict the 
behaviour of the condensate, or its excitations. Therefore, all predictions of the 
ID or 2D GPE presented in this book are necessarily only valid in appropriately 
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elongated geometries which are, however, sujfi,ciently far both from the regime of 
quasi-lD (quasi-2D) behaviour and from Tc- Finally, we note that the very rich 
ID and 2D strongly-interacting regimes can also not be described by the GPE. 
As a result, one should generally use the Gross-Pitaevskii Equation with caution 
when describing the dynamics of macroscopic excitations (solitons, vortices), for 
systems of low dimensionality at nonzero temperatures. 
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